function [FID_Grid,Action,Age,AgeState,States,pproc,DS] = prepdata_boot(opt)

% Load data

load([opt.DataFolder opt.MicroData])

% Clean nonconforming values for transportation cost:
DS = DS(double(DS(:,opt.TransCostVar)) ~= -9999,:);

% Sample selection:
% Select random subsample of fields:
if opt.SelSubSample
    [DS] = subsample(DS,opt.SubSampleFrac);
end

if opt.Buffer40km
    [DS] = DS(DS.Buffer40km == 1,:);
end

if opt.SPonly
    [DS] = DS(DS.CanasatDataset == 1,:);
end

if opt.NotPantanal
    [DS] = DS(DS.ECO_NUM ~= 7,:);
end

if opt.NotUcsIndian
    [DS] = DS(DS.indian == 0 & DS.ucs == 0,:);
end

% Eliminate fields in BA and ES:
state = floor(DS.munic*10^-5);
DS = DS(state ~= 29 & state ~= 32,:);

%% NONPARAMETRIC SAMPLE BOOTSTRAP
n       = size(DS,1);
sel     = ceil(n*unifrnd(0,1,n,1));	                                         
DS = DS(sel,:); 
DS.FID_Grid = [1:n]'; % New FID_Grid (if keep the same as before we will have repeated IDs)

%% Create action variable
% Get field ages:
[Age] = FieldAge(DS);

% Compute Action or decision variables based on observed field specific
% state change.

% Action:
Action = nan(size(Age));
AgeState = nan(size(Age));

% Actions for fields not in sugarcane:
% Keep out: (Action = 0)
Nper = size(Action,2);
for t=2:Nper
    
    % Actions for fields not in sugarcane:
    % Keep out: (Action = 0)
    Action(Age(:,t) == 0 & Age(:,t-1) == 0,t) = 0;
    % Plant sugarcane: (Action = 1)
    Action(Age(:,t) == 1 & Age(:,t-1) == 0,t) = 1;
    
    % Actions for sugarcane fields:
    % Keep field: (Action = 0)
    Action(Age(:,t) > 1 & (Age(:,t)-Age(:,t-1) == 1),t) = 0;
    % Replant field: (Action = 1)
    Action(Age(:,t) == 1 & Age(:,t-1) > 0,t) = 1;
    % Abandon sugarcane: (Action = -1)
    Action(Age(:,t) == 0 & Age(:,t-1) > 0,t) = -1;
    
    % Age state that corresponds 
    AgeState(:,t) = Age(:,t-1);
    
end

% AgeState goes from 1 to opt.MaxAge:
AgeState = AgeState + 1;

% Truncate at very old age:
if opt.TrunkAge
    AgeState( any(AgeState >= opt.TrunkAgeLevel,2), : ) = NaN;
    if opt.ElimFrac
        RandFields = rand(size(AgeState,1),1);
        ElimFrac = ((any(AgeState >= 5,2) & ~(any(Action==1,2) | any(Action==-1,2))) & (RandFields < 0.4));
        AgeState( ElimFrac, : ) = NaN;
    end
end

% Cap Age  State according to options:
AgeState(AgeState >= opt.MaxAge) = opt.MaxAge;



% Generate and classify data in states. Also run AR model for outside
% option return 
PriceChange.on = 0; % Set to off (this is ON only for counterfactual)
[States,pproc] = GenStateVariables(DS,opt,AgeState,PriceChange);



% Some fields will have no classification in terms of outside option.
% Eliminate those fields below:
Fields2keep = (States.RegionsForExoStates > 0);
disp([' Fields with no outside option classification: ']);
disp([' Total (Fraction):         '  num2str(sum(~Fields2keep)) ' (' num2str(sum(~Fields2keep)/length(Fields2keep)) ')']);
disp([' Total of S (Fraction):    '  num2str(sum(~Fields2keep & any(AgeState>1,2))) ' (' num2str(sum(~Fields2keep & any(AgeState>1,2))/sum(any(AgeState>1,2))) ')']);
disp([' '])


States.ExoTimeStates  = States.ExoTimeStates(Fields2keep,:);
States.ExoFixedStates = States.ExoFixedStates(Fields2keep,:);
States.SingleStateIdx = States.SingleStateIdx(Fields2keep,:);
States.RegionsForExoStates = States.RegionsForExoStates(Fields2keep > 0);
DS                    = DS(Fields2keep > 0,:);
AgeState              = AgeState(Fields2keep,:);
Action                = Action(Fields2keep,:);
Age                   = Age(Fields2keep,:);

AgeState              = int8(AgeState);
Action                = int8(Action);
Age                   = int8(Age);
%FID_Grid = DS.FID_Grid;
FID_Grid = uint32(DS.FID_Grid);

disp([' '])
disp([' Number of fields per outside option state (% of those in S) '])
disp('-------------------------------------------------------')
disp('  Class    OutReturn       #fields       frac. in sug. ')
disp('-------------------------------------------------------')
for r = unique(States.RegionsForExoStates)'
    AverageOut     = (pproc.outside{r}.Cte/(1 - pproc.outside{r}.AR))*1000;
    TotalInClass   = sum(States.RegionsForExoStates == r);
    TotalInClassNS = sum(States.RegionsForExoStates == r & any(AgeState>1,2));
    
    fprintf('%4d',r); fprintf('%15.2f',AverageOut); fprintf('%15d',TotalInClass); fprintf('%15.5f',TotalInClassNS/TotalInClass); fprintf('\n');
end



% % Variables from Agregate data to keep:
% VarAgregKeep = {'geocodigo', 'Year', 'micro_reg', 'meso_reg', 'uf', ... % Identifiers
%                         'p_sugarbr_lag', 'p_sugarbr', 'sb1', 'sb12', 'sb1brl', 'sb12brl', ... % Sugar price data
%                         'sb1_lag', 'sb12_lag', 'sb1brl_lag', 'sb12brl_lag', ... 
%                         'c1brl_lag', ... % Corn prices
%                         's1brl_lag', ... % Soy prices
%                         'cs_corn_qty', 'cs_corn_area', ...  % Census 2006 corn info
%                         'cs_soy_qty', 'cs_soy_area'}; % Census 2006 soy info
%                         
% Agregate.Properties.VarNames{3} = 'Year';
% Agregate.Properties.VarNames{1} = 'geocodigo';
% 
% AgregateCS = Agregate(:,VarAgregKeep);
% 
% AgregateCS = sortrows(AgregateCS,{'geocodigo','Year'});
% AgregateCS = AgregateCS(AgregateCS.Year >= 2003,:);

